% this code is for Jumps and Post-FOMC Announcement Returns in Currency Markets
% this code is making a figure showing return dynamics depending on jump
% and direction of monetary policies and surprises (figures 2 & 5)

% figure 4 can be drwan by MS Excel with raw data

clear
clc

load(using 15 min filtered return data)

% 1) figure 2

fomcrxeease=nan(npday*4,nc+1,160);
fomcrxetight=nan(npday*4,nc+1,160);
fomcrxeneutral=nan(npday*4,nc+1,160);

tpone=0;

for j=1:nt-1
    if usfomc(j,1)==1 && usfomc(j+1,1)==1       % fomc indicator
        tpone=tpone+1;
        if usfomc(j,2)==1   %ease?
            fomcrxeease(:,:,tpone)=fomcrxe(j-npday*2+1:j+npday*2,:);
        end
        if usfomc(j,3)==1   %neutral?
            fomcrxeneutral(:,:,tpone)=fomcrxe(j-npday*2+1:j+npday*2,:);
        end
        if usfomc(j,4)==1   %tight?
            fomcrxetight(:,:,tpone)=fomcrxe(j-npday*2+1:j+npday*2,:);
        end
    end
end
clear j tpone

% indicator matrix for negative jump

jdncub=nan(npday*4,nc+1,160);
jdnnocub=nan(npday*4,nc+1,160);

tpone=0;

for j=1:nt-1
    if usfomc(j,1)==1 && usfomc(j+1,1)==1
        tpone=tpone+1;
        tpjp=jdn(j+1:j+npday/2,:);
        tpsumjp=nansum(tpjp,1);
        for i=1:nc+1
            if tpsumjp(1,i)>0
                jdncub(:,i,tpone)=1;
            end
            if tpsumjp(1,i)<1
                jdnnocub(:,i,tpone)=1;
            end
        end
        clear i
    end
end
clear j tpone

fomcrxeeasenjp=nan(npday*4,nc+1,160);
fomcrxetightnjp=nan(npday*4,nc+1,160);
fomcrxeneutralnjp=nan(npday*4,nc+1,160);

fomcrxeeasenonjp=nan(npday*4,nc+1,160);
fomcrxetightnonjp=nan(npday*4,nc+1,160);
fomcrxeneutralnonjp=nan(npday*4,nc+1,160);

for k=1:160
    fomcrxeeasenjp(:,:,k)=fomcrxeease(:,:,k).*jdncub(:,:,k);
    fomcrxetightnjp(:,:,k)=fomcrxetight(:,:,k).*jdncub(:,:,k);
    fomcrxeneutralnjp(:,:,k)=fomcrxeneutral(:,:,k).*jdncub(:,:,k);
    
    fomcrxeeasenonjp(:,:,k)=fomcrxeease(:,:,k).*jdnnocub(:,:,k);
    fomcrxetightnonjp(:,:,k)=fomcrxetight(:,:,k).*jdnnocub(:,:,k);
    fomcrxeneutralnonjp(:,:,k)=fomcrxeneutral(:,:,k).*jdnnocub(:,:,k);
end
clear k

sumfomcrxeeasenjp=nan(npday*4,2,160);
sumfomcrxetightnjp=nan(npday*4,2,160);
sumfomcrxeneutralnjp=nan(npday*4,2,160);

sumfomcrxeeasenonjp=nan(npday*4,2,160);
sumfomcrxetightnonjp=nan(npday*4,2,160);
sumfomcrxeneutralnonjp=nan(npday*4,2,160);

for k=1:160
    sumfomcrxeeasenjp(:,:,k)=[nanmean(fomcrxeeasenjp(:,1:nc,k),2) ...
        fomcrxeeasenjp(:,nc+1,k)];
    sumfomcrxetightnjp(:,:,k)=[nanmean(fomcrxetightnjp(:,1:nc,k),2) ...
        fomcrxetightnjp(:,nc+1,k)];
    sumfomcrxeneutralnjp(:,:,k)=[nanmean(fomcrxeneutralnjp(:,1:nc,k),2) ...
        fomcrxeneutralnjp(:,nc+1,k)];

    sumfomcrxeeasenonjp(:,:,k)=[nanmean(fomcrxeeasenonjp(:,1:nc,k),2) ...
        fomcrxeeasenonjp(:,nc+1,k)];
    sumfomcrxetightnonjp(:,:,k)=[nanmean(fomcrxetightnonjp(:,1:nc,k),2) ...
        fomcrxetightnonjp(:,nc+1,k)];
    sumfomcrxeneutralnonjp(:,:,k)=[nanmean(fomcrxeneutralnonjp(:,1:nc,k),2) ...
        fomcrxeneutralnonjp(:,nc+1,k)];
end
clear k

repfomcrxeeasenjp=nan(npday*4,2);
repfomcrxetightnjp=nan(npday*4,2);
repfomcrxeneutralnjp=nan(npday*4,2);

repfomcrxeeasenonjp=nan(npday*4,2);
repfomcrxetightnonjp=nan(npday*4,2);
repfomcrxeneutralnonjp=nan(npday*4,2);

for j=1:npday*4
    repfomcrxeeasenjp(j,1)=nanmean(sumfomcrxeeasenjp(j,1,:));
    repfomcrxeeasenjp(j,2)=nanmean(sumfomcrxeeasenjp(j,2,:));
    
    repfomcrxetightnjp(j,1)=nanmean(sumfomcrxetightnjp(j,1,:));
    repfomcrxetightnjp(j,2)=nanmean(sumfomcrxetightnjp(j,2,:));
    
    repfomcrxeneutralnjp(j,1)=nanmean(sumfomcrxeneutralnjp(j,1,:));
    repfomcrxeneutralnjp(j,2)=nanmean(sumfomcrxeneutralnjp(j,2,:));
    
    repfomcrxeeasenonjp(j,1)=nanmean(sumfomcrxeeasenonjp(j,1,:));
    repfomcrxeeasenonjp(j,2)=nanmean(sumfomcrxeeasenonjp(j,2,:));
    
    repfomcrxetightnonjp(j,1)=nanmean(sumfomcrxetightnonjp(j,1,:));
    repfomcrxetightnonjp(j,2)=nanmean(sumfomcrxetightnonjp(j,2,:));
    
    repfomcrxeneutralnonjp(j,1)=nanmean(sumfomcrxeneutralnonjp(j,1,:));
    repfomcrxeneutralnonjp(j,2)=nanmean(sumfomcrxeneutralnonjp(j,2,:));
end
clear j

repfomccumrxeeasenjp=zeros(npday*4+1,2);
repfomccumrxetightnjp=zeros(npday*4+1,2);
repfomccumrxeneutralnjp=zeros(npday*4+1,2);

repfomccumrxeeasenonjp=zeros(npday*4+1,2);
repfomccumrxetightnonjp=zeros(npday*4+1,2);
repfomccumrxeneutralnonjp=zeros(npday*4+1,2);

for j=1:npday*4
    repfomccumrxeeasenjp(j+1,:)=repfomccumrxeeasenjp(j,:)...
        +repfomcrxeeasenjp(j,:);
    repfomccumrxetightnjp(j+1,:)=repfomccumrxetightnjp(j,:)...
        +repfomcrxetightnjp(j,:);
    repfomccumrxeneutralnjp(j+1,:)=repfomccumrxeneutralnjp(j,:)...
        +repfomcrxeneutralnjp(j,:);
    
    repfomccumrxeeasenonjp(j+1,:)=repfomccumrxeeasenonjp(j,:)...
        +repfomcrxeeasenonjp(j,:);
    repfomccumrxetightnonjp(j+1,:)=repfomccumrxetightnonjp(j,:)...
        +repfomcrxetightnonjp(j,:);
    repfomccumrxeneutralnonjp(j+1,:)=repfomccumrxeneutralnonjp(j,:)...
        +repfomcrxeneutralnonjp(j,:);
end
clear j

repfomccumrxecomp=[repfomccumrxeeasenjp(:,1) repfomccumrxeeasenonjp(:,1) ...
    repfomccumrxetightnjp(:,1) repfomccumrxetightnonjp(:,1)];

% 2) figure 5

% 2.1) making return 3d matrix

temp=0;
fomcrxecub=nan(npday*4,nc+1,160);
for j=1:nt
    if fomcwindow(j,1)==1 && (fomcwindow(j-1,1)>=0)==0
        temp=temp+1;
        fomcrxecub(:,:,temp)=fomcrxe(j:j+npday*4-1,:);
    end
end
clear i temp

% 2.2) classifying posi>zero>zega policy surprises

surpp=nan(npday*4,nc+1,160);
surpp(:,:,4)=ones(npday*4,nc+1);
surpp(:,:,5)=ones(npday*4,nc+1);
surpp(:,:,7)=ones(npday*4,nc+1);
surpp(:,:,8)=ones(npday*4,nc+1);
surpp(:,:,10)=ones(npday*4,nc+1);
surpp(:,:,11)=ones(npday*4,nc+1);
surpp(:,:,44)=ones(npday*4,nc+1);
surpp(:,:,45)=ones(npday*4,nc+1);
surpp(:,:,46)=ones(npday*4,nc+1);
surpp(:,:,47)=ones(npday*4,nc+1);
surpp(:,:,48)=ones(npday*4,nc+1);
surpp(:,:,50)=ones(npday*4,nc+1);
surpp(:,:,52)=ones(npday*4,nc+1);
surpp(:,:,54)=ones(npday*4,nc+1);
surpp(:,:,55)=ones(npday*4,nc+1);
surpp(:,:,56)=ones(npday*4,nc+1);
surpp(:,:,57)=ones(npday*4,nc+1);
surpp(:,:,58)=ones(npday*4,nc+1);
surpp(:,:,60)=ones(npday*4,nc+1);
surpp(:,:,78)=ones(npday*4,nc+1);
surpp(:,:,80)=ones(npday*4,nc+1);
surpp(:,:,81)=ones(npday*4,nc+1);
surpp(:,:,83)=ones(npday*4,nc+1);
surpp(:,:,85)=ones(npday*4,nc+1);
surpp(:,:,86)=ones(npday*4,nc+1);
surpp(:,:,87)=ones(npday*4,nc+1);
surpp(:,:,88)=ones(npday*4,nc+1);
surpp(:,:,89)=ones(npday*4,nc+1);
surpp(:,:,90)=ones(npday*4,nc+1);
surpp(:,:,92)=ones(npday*4,nc+1);
surpp(:,:,93)=ones(npday*4,nc+1);
surpp(:,:,94)=ones(npday*4,nc+1);
surpp(:,:,95)=ones(npday*4,nc+1);
surpp(:,:,96)=ones(npday*4,nc+1);
surpp(:,:,97)=ones(npday*4,nc+1);
surpp(:,:,98)=ones(npday*4,nc+1);
surpp(:,:,99)=ones(npday*4,nc+1);
surpp(:,:,100)=ones(npday*4,nc+1);
surpp(:,:,101)=ones(npday*4,nc+1);
surpp(:,:,102)=ones(npday*4,nc+1);
surpp(:,:,103)=ones(npday*4,nc+1);
surpp(:,:,104)=ones(npday*4,nc+1);
surpp(:,:,105)=ones(npday*4,nc+1);
surpp(:,:,106)=ones(npday*4,nc+1);
surpp(:,:,107)=ones(npday*4,nc+1);
surpp(:,:,108)=ones(npday*4,nc+1);
surpp(:,:,109)=ones(npday*4,nc+1);
surpp(:,:,110)=ones(npday*4,nc+1);
surpp(:,:,111)=ones(npday*4,nc+1);
surpp(:,:,112)=ones(npday*4,nc+1);
surpp(:,:,113)=ones(npday*4,nc+1);
surpp(:,:,114)=ones(npday*4,nc+1);
surpp(:,:,115)=ones(npday*4,nc+1);
surpp(:,:,116)=ones(npday*4,nc+1);
surpp(:,:,117)=ones(npday*4,nc+1);
surpp(:,:,118)=ones(npday*4,nc+1);
surpp(:,:,119)=ones(npday*4,nc+1);
surpp(:,:,120)=ones(npday*4,nc+1);
surpp(:,:,121)=ones(npday*4,nc+1);
surpp(:,:,122)=ones(npday*4,nc+1);
surpp(:,:,123)=ones(npday*4,nc+1);
surpp(:,:,124)=ones(npday*4,nc+1);
surpp(:,:,125)=ones(npday*4,nc+1);
surpp(:,:,126)=ones(npday*4,nc+1);
surpp(:,:,127)=ones(npday*4,nc+1);
surpp(:,:,128)=ones(npday*4,nc+1);
surpp(:,:,129)=ones(npday*4,nc+1);
surpp(:,:,130)=ones(npday*4,nc+1);
surpp(:,:,131)=ones(npday*4,nc+1);
surpp(:,:,132)=ones(npday*4,nc+1);
surpp(:,:,133)=ones(npday*4,nc+1);
surpp(:,:,134)=ones(npday*4,nc+1);
surpp(:,:,135)=ones(npday*4,nc+1);
surpp(:,:,136)=ones(npday*4,nc+1);
surpp(:,:,137)=ones(npday*4,nc+1);
surpp(:,:,138)=ones(npday*4,nc+1);
surpp(:,:,139)=ones(npday*4,nc+1);
surpp(:,:,140)=ones(npday*4,nc+1);
surpp(:,:,141)=ones(npday*4,nc+1);
surpp(:,:,142)=ones(npday*4,nc+1);
surpp(:,:,143)=ones(npday*4,nc+1);

surpz=nan(npday*4,nc+1,160);
surpz(:,:,1)=ones(npday*4,nc+1);
surpz(:,:,3)=ones(npday*4,nc+1);
surpz(:,:,6)=ones(npday*4,nc+1);
surpz(:,:,9)=ones(npday*4,nc+1);
surpz(:,:,12)=ones(npday*4,nc+1);
surpz(:,:,13)=ones(npday*4,nc+1);
surpz(:,:,14)=ones(npday*4,nc+1);
surpz(:,:,15)=ones(npday*4,nc+1);
surpz(:,:,16)=ones(npday*4,nc+1);
surpz(:,:,22)=ones(npday*4,nc+1);
surpz(:,:,25)=ones(npday*4,nc+1);
surpz(:,:,26)=ones(npday*4,nc+1);
surpz(:,:,27)=ones(npday*4,nc+1);
surpz(:,:,28)=ones(npday*4,nc+1);
surpz(:,:,29)=ones(npday*4,nc+1);
surpz(:,:,30)=ones(npday*4,nc+1);
surpz(:,:,32)=ones(npday*4,nc+1);
surpz(:,:,33)=ones(npday*4,nc+1);
surpz(:,:,34)=ones(npday*4,nc+1);
surpz(:,:,35)=ones(npday*4,nc+1);
surpz(:,:,37)=ones(npday*4,nc+1);
surpz(:,:,38)=ones(npday*4,nc+1);
surpz(:,:,39)=ones(npday*4,nc+1);
surpz(:,:,40)=ones(npday*4,nc+1);
surpz(:,:,41)=ones(npday*4,nc+1);
surpz(:,:,42)=ones(npday*4,nc+1);
surpz(:,:,43)=ones(npday*4,nc+1);
surpz(:,:,49)=ones(npday*4,nc+1);
surpz(:,:,51)=ones(npday*4,nc+1);
surpz(:,:,53)=ones(npday*4,nc+1);
surpz(:,:,59)=ones(npday*4,nc+1);
surpz(:,:,61)=ones(npday*4,nc+1);
surpz(:,:,62)=ones(npday*4,nc+1);
surpz(:,:,63)=ones(npday*4,nc+1);
surpz(:,:,64)=ones(npday*4,nc+1);
surpz(:,:,65)=ones(npday*4,nc+1);
surpz(:,:,66)=ones(npday*4,nc+1);
surpz(:,:,67)=ones(npday*4,nc+1);
surpz(:,:,68)=ones(npday*4,nc+1);
surpz(:,:,69)=ones(npday*4,nc+1);
surpz(:,:,72)=ones(npday*4,nc+1);
surpz(:,:,76)=ones(npday*4,nc+1);
surpz(:,:,77)=ones(npday*4,nc+1);
surpz(:,:,79)=ones(npday*4,nc+1);
surpz(:,:,82)=ones(npday*4,nc+1);
surpz(:,:,84)=ones(npday*4,nc+1);
surpz(:,:,91)=ones(npday*4,nc+1);
surpz(:,:,144)=ones(npday*4,nc+1);
surpz(:,:,146)=ones(npday*4,nc+1);
surpz(:,:,148)=ones(npday*4,nc+1);
surpz(:,:,152)=ones(npday*4,nc+1);
surpz(:,:,154)=ones(npday*4,nc+1);
surpz(:,:,158)=ones(npday*4,nc+1);
surpz(:,:,160)=ones(npday*4,nc+1);

surpn=nan(npday*4,nc+1,160);
surpn(:,:,2)=ones(npday*4,nc+1);
surpn(:,:,17)=ones(npday*4,nc+1);
surpn(:,:,18)=ones(npday*4,nc+1);
surpn(:,:,19)=ones(npday*4,nc+1);
surpn(:,:,20)=ones(npday*4,nc+1);
surpn(:,:,21)=ones(npday*4,nc+1);
surpn(:,:,23)=ones(npday*4,nc+1);
surpn(:,:,24)=ones(npday*4,nc+1);
surpn(:,:,31)=ones(npday*4,nc+1);
surpn(:,:,36)=ones(npday*4,nc+1);
surpn(:,:,70)=ones(npday*4,nc+1);
surpn(:,:,71)=ones(npday*4,nc+1);
surpn(:,:,73)=ones(npday*4,nc+1);
surpn(:,:,74)=ones(npday*4,nc+1);
surpn(:,:,75)=ones(npday*4,nc+1);
surpn(:,:,145)=ones(npday*4,nc+1);
surpn(:,:,147)=ones(npday*4,nc+1);
surpn(:,:,149)=ones(npday*4,nc+1);
surpn(:,:,150)=ones(npday*4,nc+1);
surpn(:,:,151)=ones(npday*4,nc+1);
surpn(:,:,153)=ones(npday*4,nc+1);
surpn(:,:,155)=ones(npday*4,nc+1);
surpn(:,:,156)=ones(npday*4,nc+1);
surpn(:,:,157)=ones(npday*4,nc+1);
surpn(:,:,159)=ones(npday*4,nc+1);

% 2.3) indicator matrix for negative jump

jdncub=nan(npday*4,nc+1,160);
jdnnocub=nan(npday*4,nc+1,160);

for k=1:160
    for i=1:nc+1
        if jdnvhlfdaywindow(k,i)>0
            jdncub(:,i,k)=ones(npday*4,1);
        end
        if jdnvhlfdaywindow(k,i)==0
            jdnnocub(:,i,k)=ones(npday*4,1);
        end
        if (jdnvhlfdaywindow(k,i)>=0)+(jdnvhlfdaywindow(k,i)<0)==0
            jdnnocub(:,i,k)=ones(npday*4,1);
        end
    end
    clear i
end
clear k

sumdetretsurppjdn=nan(npday*4,nc+1,160);
sumdetretsurppjdnno=nan(npday*4,nc+1,160);
sumdetretsurpzjdn=nan(npday*4,nc+1,160);
sumdetretsurpzjdnno=nan(npday*4,nc+1,160);
sumdetretsurpnjdn=nan(npday*4,nc+1,160);
sumdetretsurpnjdnno=nan(npday*4,nc+1,160);

sumdetavgretsurppjdn=nan(npday*4,nc+1);
sumdetavgretsurppjdnno=nan(npday*4,nc+1);
sumdetavgretsurpzjdn=nan(npday*4,nc+1);
sumdetavgretsurpzjdnno=nan(npday*4,nc+1);
sumdetavgretsurpnjdn=nan(npday*4,nc+1);
sumdetavgretsurpnjdnno=nan(npday*4,nc+1);

for k=1:160
    sumdetretsurppjdn(:,:,k)=fomcrxecub(:,:,k).*surpp(:,:,k).*jdncub(:,:,k);
    sumdetretsurppjdnno(:,:,k)=fomcrxecub(:,:,k).*surpp(:,:,k).*jdnnocub(:,:,k);
    sumdetretsurpzjdn(:,:,k)=fomcrxecub(:,:,k).*surpz(:,:,k).*jdncub(:,:,k);
    sumdetretsurpzjdnno(:,:,k)=fomcrxecub(:,:,k).*surpz(:,:,k).*jdnnocub(:,:,k);
    sumdetretsurpnjdn(:,:,k)=fomcrxecub(:,:,k).*surpn(:,:,k).*jdncub(:,:,k);
    sumdetretsurpnjdnno(:,:,k)=fomcrxecub(:,:,k).*surpn(:,:,k).*jdnnocub(:,:,k);
end
clear k

for i=1:nc+1
    for j=1:npday*4
        tpone=nan(160,1);
        tptwo=nan(160,1);
        tpthr=nan(160,1);
        tpfor=nan(160,1);
        tpfif=nan(160,1);
        tpsix=nan(160,1);
        for k=1:160
            tpone(k,1)=sumdetretsurppjdn(j,i,k);
            tptwo(k,1)=sumdetretsurppjdnno(j,i,k);
            tpthr(k,1)=sumdetretsurpzjdn(j,i,k);
            tpfor(k,1)=sumdetretsurpzjdnno(j,i,k);
            tpfif(k,1)=sumdetretsurpnjdn(j,i,k);
            tpsix(k,1)=sumdetretsurpnjdnno(j,i,k);
        end
        clear k
        sumdetavgretsurppjdn(j,i)=nanmean(tpone,1);
        sumdetavgretsurppjdnno(j,i)=nanmean(tptwo,1);
        sumdetavgretsurpzjdn(j,i)=nanmean(tpthr,1);
        sumdetavgretsurpzjdnno(j,i)=nanmean(tpfor,1);
        sumdetavgretsurpnjdn(j,i)=nanmean(tpfif,1);
        sumdetavgretsurpnjdnno(j,i)=nanmean(tpsix,1);
    end
    clear j
end
clear i

tpone=sumdetavgretsurppjdn;
tptwo=sumdetavgretsurppjdnno;
tpthr=sumdetavgretsurpzjdn;
tpfor=sumdetavgretsurpzjdnno;
tpfif=sumdetavgretsurpnjdn;
tpsix=sumdetavgretsurpnjdnno;

for i=1:nc+1
    for j=1:npday*4
        if (sumdetavgretsurppjdn(j,i)>=0)+(sumdetavgretsurppjdn(j,i)<0)==0
            tpone(j,i)=0;
        end
        if (sumdetavgretsurppjdnno(j,i)>=0)+(sumdetavgretsurppjdnno(j,i)<0)==0
            tptwo(j,i)=0;
        end
        if (sumdetavgretsurpzjdn(j,i)>=0)+(sumdetavgretsurpzjdn(j,i)<0)==0
            tpthr(j,i)=0;
        end
        if (sumdetavgretsurpzjdnno(j,i)>=0)+(sumdetavgretsurpzjdnno(j,i)<0)==0
            tpfor(j,i)=0;
        end
        if (sumdetavgretsurpnjdn(j,i)>=0)+(sumdetavgretsurpnjdn(j,i)<0)==0
            tpfif(j,i)=0;
        end
        if (sumdetavgretsurpnjdnno(j,i)>=0)+(sumdetavgretsurpnjdnno(j,i)<0)==0
            tpsix(j,i)=0;
        end
    end
    clear j
end
clear i

sumdetcumretsurppjdn=zeros(npday*4+1,nc+1);
sumdetcumretsurppjdnno=zeros(npday*4+1,nc+1);
sumdetcumretsurpzjdn=zeros(npday*4+1,nc+1);
sumdetcumretsurpzjdnno=zeros(npday*4+1,nc+1);
sumdetcumretsurpnjdn=zeros(npday*4+1,nc+1);
sumdetcumretsurpnjdnno=zeros(npday*4+1,nc+1);

for j=2:npday*4+1
    sumdetcumretsurppjdn(j,:)=sumdetcumretsurppjdn(j-1,:)+tpone(j-1,:);
    sumdetcumretsurppjdnno(j,:)=sumdetcumretsurppjdnno(j-1,:)+tptwo(j-1,:);
    sumdetcumretsurpzjdn(j,:)=sumdetcumretsurpzjdn(j-1,:)+tpthr(j-1,:);
    sumdetcumretsurpzjdnno(j,:)=sumdetcumretsurpzjdnno(j-1,:)+tpfor(j-1,:);
    sumdetcumretsurpnjdn(j,:)=sumdetcumretsurpnjdn(j-1,:)+tpfif(j-1,:);
    sumdetcumretsurpnjdnno(j,:)=sumdetcumretsurpnjdnno(j-1,:)+tpsix(j-1,:);
end
clear j

sumdetcumretsurppjdn(:,2)=[];
sumdetcumretsurppjdnno(:,2)=[];
sumdetcumretsurpzjdn(:,2)=[];
sumdetcumretsurpzjdnno(:,2)=[];
sumdetcumretsurpnjdn(:,2)=[];
sumdetcumretsurpnjdnno(:,2)=[];

% for report
aggdetcumretsurppjdn=nanmean(sumdetcumretsurppjdn(:,1:nc-1),2);
aggdetcumretsurppjdnno=nanmean(sumdetcumretsurppjdnno(:,1:nc-1),2);
aggdetcumretsurpzjdn=nanmean(sumdetcumretsurpzjdn(:,1:nc-1),2);
aggdetcumretsurpzjdnno=nanmean(sumdetcumretsurpzjdnno(:,1:nc-1),2);
aggdetcumretsurpnjdn=nanmean(sumdetcumretsurpnjdn(:,1:nc-1),2);
aggdetcumretsurpnjdnno=nanmean(sumdetcumretsurpnjdnno(:,1:nc-1),2);
